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Abstract. We study the work statistics of a periodically-driven integrable closed 
quantum system, addressing in particular the role played by the presence of 
a quantum critical point. Taking the example of a one-dimensional transverse 
Ising model in the presence of a spatially homogeneous but periodically time- 
varying transverse field of frequency ljq, we arrive at the characteristic cumulant 
generating function G(u), which is then used to calculate the work distribution 
function P{W). By applying the Floquet theory we show that, in the infinite time 
limit, P{W) converges, at zero temperature, towards an asymptotic steady state 
value whose small-IV behaviour depends only on the properties of the small-wave- 
vector modes and on a few important ingredients: the time-averaged value of the 
transverse field, ho, the initial transverse field, hj, and the equilibrium quantum 
critical point he, which we find to generate a sequence of non-equilibrium critical 
points = he + Iloo/2^ with I integer. When hi 7^ he, we find a “universal” edge 
singularity in P{W) at a threshold value of Wth = 2|hi — he\ which is entirely 
determined by hi. The form of that singularity — Dirac delta derivative or square 
root — depends on ho being or not at a non-equilibrium critical point h*^. On 
the contrary, when hi = he, G{u) decays as a power-law for large u, leading to 
different types of edge singularity at Wth = 0. Generalizing our calculations to 
the finite temperature case, the irreversible entropy generated by the periodic 
driving is also shown to reach a steady state value in the infinite time limit. 
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1. Introduction 

In recent years, there have been many theoretical studies aimed at understanding 
the non-equilibrium dynamics of closed quantum systems [U [5], inspired by a series of 
experiments on cold atomic gases, which are nearly isolated systems with long phase- 
coherence times, allowing for the study of a coherent quantum dynamics over long 
time scales 011 IB El- These experiments have paved the way to addressing many 
fundamental questions, such as the role of integrability in thermalization following a 
quench [T] or the universal scaling of the defects generated when a system is driven 
across a quantum critical point 0[7]. Usually, a quantum system is driven out of 
equilibrium by a slow ramping (annealing) or by a sudden quench of a parameter 
of the Hamiltonian (for example, the transverse magnetic field in the quantum Ising 
model discussed in this paper); the subsequent dynamical response of the system is 
encoded in several quantities, e.g., the Loschmidt echo 0ini[ini[iii[n], dynamical 
correlation functions m, the growth of entanglement entropy m, or the change 
in diagonal entropy m- In parallel, in the context of the Jarzynski equality [I6| 
and non-equilibrium fluctuation relations m, the question of the emergence of 
thermodynamical laws from a finite quantum system driven out of equilibrium and the 
generation of irreversible entropy have been addressed in several recent works [HUS]- 
One of the ways to characterize the dynamics of an out-of-equilibrium quantum 
system is to explore the statistics of the performed work [20l [211 Il2] • Given the non¬ 
equilibrium nature of the driving protocol, the work (W) is a stochastic variable and 
hence described by a probability distribution P{W). Recently, from the analysis of the 
work statistics in systems quenched across a critical point, it has been shown that non- 
analyiticities in P(W = 0) at some critical times can mark the existence of a sequence 
of dynamical phase transitions in real time [23] • Furthermore, the knowledge of P{W) 
enables us to obtain information about some universal features, by connecting it to 
the critical Casimir effect, for sudden quenches ending near the critical point |21j : in 


particular, there exists a power-law edge singularity in P{W), for small W, which is 
characterized by an exponent that is independent of the choice of protocol, but rather 
depends just on the initial and final values of the parameter being quenched |24j . It is 
usually convenient to define and work with the characteristic function G{u), obtained 
by Fourier transforming P{W). The characteristic function G{u) has been shown to 
be closely related to the Loschmidt echo of a quenched quantum system, both at zero 
[a ED] and finite temperature [25] . 

In this work we focus on a periodically driven integrable closed quantum system, 
namely the transverse Ising chain EiET], with a spatially homogeneous but time- 
periodic transverse field h{t) = h{t + 2Tr/oJo), and calculate P{W) stroboscopically at 
the end of n complete periods r = 2ttIujq. It has been shown in the literature [28l ED] 
that the system reaches a periodic steady state in the limit n —>■ oo: the residual 
energy (which is in fact the first moment of P{W)) — and indeed essentially any 
local observable — reaches a stationary value, when observed at times = nr, in the 
thermodynamic limit. We have observed a similar relaxation to a periodic steady 
condition also for a genuinely quantum non-local object, the so-called dynamical 
hdelity |30| . In the present work we will investigate what happens to G(u), and 
hence to P{W), under such periodic driving in the asymptotic limit n —?> oo, and 
address the question of the universal behaviour emerging in the small-VF region. 

Working within the framework of the Floquet theory iniisi], and assuming that 
the initial state is a Gibbs state at temperature T, we provide an analytical form of 
the characteristic function G{u) in terms of Floquet quasi-energies and corresponding 
overlaps between the initial state and the Floquet eigenstates. G{u) is then used 
to arrive at an exact expression of P(W) at zero temperature: we demonstrate 
that indeed P{W) also tends to “synchronize” with the periodic driving in the 
limit n —)■ 00 , converging to a well-defined asymptotic work distribution function 
Poo{W), whose small-kF behaviour depends only on the properties of the small- 
wavevector modes, ultimately controlled by the time-averaged value of the transverse 
held Hq = T~^ fj h(t)dt, and its initial value = h{t = 0). All universal features of 
the work distribution are encoded into a singularity of P{W) for small W. The position 
Wth = 21 hi — he I of the edge of that singularity depends only on the distance of hi from 
the equilibrium critical point he, while its detailed form — a Dirac delta derivative or 
a square root singularity — depends on hg being or not at a non-equilibrium critical 
point h^i — determined by h^,i = he + Iojq/2 with I integer — where the Floquet 
spectrum turns out to be gapless. When = he, the generating function G{u) has a 
power law decay l/{iu)P, which results into a step singularity in Poo(bF) at IFth = 0 
for p = I, or a mild increase for larger integer p. The non-equilibrium phase 

transitions we hnd are reminescent of those discussed in Ref. |33| , but are different in 
many aspects, notably in residing at low frequency, as we will discuss. 

An important aspect of these results is that the small-IF properties of Poo{W) are 
determined (as we are going to show in detail in the paper) exclusively by the small-fc 
modes: the aspects of the dynamics independent of the details of the driving protocol 
rely only on the large wavelenght modes which encode the universal properties of the 
ground state in the static system |34l |35| . Remarkably, a similar fact can be observed 
in two coupled Luttinger liquids undergoing a quantum quench [3^: if the quenched 
operator is relevant in the renormalization group sense (and then affecting mainly 
the long wavelength modes) then the phase coherence evolves with a universal scaling 
function. 

Finally we move to the finite temperature case: starting from an initial mixed 


state at finite temperature, we will show that the irreversible entropy generated during 
the periodic driving tends to “synchronize” with the driving for n —> oo. It saturates 
to a steady state value for large wq, displaying a sequence of dips and peaks for small 
and intermediate values of uiq, respectively. We note in passing a recent study on 
the work statistics of a periodically driven system which has explored the universal 
properties of the rate function, that is found to satisfy a lower bound and has a zero 
when W matches the residual energy m- 

The structure of the paper is as follows. In Sec. we summarize the basic 
definitions and properties of P{W) and G{u). Next, in Sec. [^we focus on the case 
of a quantum Ising chain undergoing a uniform generic periodic driving, showing 
in Sec, [i] that the stroboscopic G{u) and P{W) converge to an asymptotic value 


(see 


Appendix A for mathematical details). In Sec.j^we discuss the behaviour of the 


asymptotic P{W) in the small-W regime, showing that it is independent of the details 
of the driving protocol; we substantiate our analysis with numerical results obtained 
in the case of a sinusoidal driving. In Sec. we discuss briefly the case with finite 
temperature T in the context of irreversible entropy generation, and in Sec. [7] we draw 
our conclusions together with perspectives of future work. Technical details of our 
derivations are contained in three appendices. 


2. The work distribution P{W) and its characteristic function G{u) 

We present here, for the readers’ convenience, some basic facts about the work 
distribution function, following Refs. [ni [Ml I3S]- Suppose that a closed quantum 
system undergoes a time-dependent driving such that its Hamiltonian is H(t), while 
the system was at the initial time to = 0 in a given (possibly mixed) state p(0). If p„(0) 
is the probability that the system has energy A„(0) at the initial time, and P{mt(\nO) 
the conditional probability that the system is observed to have energy Em{t{) at some 
later time tf, then the work distribution function is defined as: 


PtAW) = ^ <5(IT - + A„(0)) P(mtf|nO) p„(0) . 


( 1 ) 


To deal with the Dirac’s deltas appearing in the definition of Ptf{W) it is convenient 
to study the Fourier transform of Pu (IT), arriving at the characteristic function: 




GtAu) = 


dWe^^^Pt,{W) . 


( 2 ) 


With simple manipulations, and introducing the unitary evolution operator U{ti,0) 
for the closed quantum system, it results that: 

Gt,(w) = Tr (l7l’(<f,0)e“^(‘f)l7(tf,0)e-™^(°)p(0)) = Tr , 

where Elnitf) = C^'*'(tf, 0)i?(tf)[/(tf, 0) is the final Hamiltonian in Heisenberg 
representation, and we have assumed that the initial state is such that [p(0), iJ(0)] = 0 
(a Gibbs or micro-canonical state would do that). We note, in passing, that the 
quantum Jarzynski equality mns] follows immediately from the previous expression 
by taking u = ij3 and assuming an initial Gibbs state p(0) = 

/ + 00 

PtAW) = . 

-OO 


( 4 ) 



3. The uniform periodically driven quantum Ising chain 


Let us now specialize our discussion to the quantum Ising chain in a uniform time- 
periodic transverse field (although some progress might be done in the general non- 
homogeneous case). The Hamiltonian we consider is |M] : 

^ T 

= -f /i(t)(Tj] , (5) 


where h{t) = h(t+T) is a generic uniform-in-space transverse field which is periodically 
driven with frequency ujq = ^tt/t around the average value 

ho = - [ h(t) dt . (6) 

T Jo 


The Hamiltonian can always be recast in the form of a quadratic fermionic model, 
thanks to the Jordan-Wigner transformation |40j . Upon Fourier transforming in space 
to the relevant Jordan-Wigner fermions, we get: 


ABC 


ABC 


H+{t) = ^ Hk{t) = ^ [ 


^-k 


IHIfc(t) 



where the 2x2 matrix ]HIfc(t) has the form: 




ek{t) -lAk 
iAk -ek{t) \ ’ 


(7) 


( 8 ) 


with efc(t) = h{t) — cosk and = sinfc, having set J = 1. Notice that the previous 
Hamiltonian is really only a part of the total H = H'^ + H~, i.e., the part living in 
the subspace with even fermion-parity, for which anti-periodic boundary conditions 
(ABC) apply, and k = {2n — l)7r/L with n = 1 ■ ■ ■ L/2 [JD]. This is certainly enough 
for describing the ground state and the dynamics starting from the ground state. At 
finite temperature, the contributions due to the extra odd-fermion-parity term, H~, 
corresponding to periodic boundary conditions fc-values, are automatically accounted 
for when we transform the sum over k into an integral over the whole Brillouin zone 
k G [0, tt]. 

We assume that the system is, at time t = 0, in the Gibbs ensemble at temperature 
T for the Hamiltonian H{0). Following a standard procedure [23113 SO]! the initial 
problem is diagonalized by introducing new fermionic operators = Ukcl + VkC_f., 
in terms of which the BCS-ground state for each k reads |^f’’) = [uk + UfeCfccLfe]|0), 
with energy -Ek = -\/e^(0) + A^, and the excited state is = 

[vk + WfcC^c^j,]|0), with energy +Ek. If we are interested in calculating Gti=nTiu) at 
integer multiples of the period t = 27r/a;o, it is sufficient to construct the Floquet 
modes (0)) with corresponding Floquet quasi-energies EfJ-k] this must be done, in 
general, by a numerical integration of the Schrodinger equation, for each k, over a 
single period |28l |30| . The outcome of that calculation provides the relevant overlaps 
= (<^J(0)IV'f), with |r+|2 -h |r-|2 = 1, and, by unitarity, {(t)k = Tir^)*■ 

With these basic ingredients, the derivation of a general expression of Gnriu) for 
the uniformly driven Ising chain follows essentially the steps outlined in Ref. m, 










generalized to an arbitrary finite temperature. The final result can be cast in the 
form: 


ABC 


In GnT{u) = In |l - ^ sin^ (/Tfcnr) (1 - e^“^'')(l - /fc) + (1 - e | , 


fc >0 


1 + 


(9) 

where fk = l/(e^^^'= +1) denotes the Fermi occupation function (observe that creating 
an excitation costs here an energy 2Ek) and 


qk = 




( 10 ) 


Had we chosen to work with the Laplace transform of the P{W), rather than with the 
Fourier transform G{u), we would have obtained a similar expression with the formal 
replacement u ^ is. For T = 0 we would get 


ABC „ 

lnG^r°(*s) = 51 Y-^sin^(^fcnT)(l-e"^'’'®")| , (11) 


which actually shows that the expression for G{is) is better behaved at T = 0, at 
the expense of having to perform an inverse Laplace transform to recover P{W). In 
the thermodynamic limit, transforming the sum over k into an integral, we can finally 
write: 


lnGL-°(^g) 

L 


^dk , 

— In ■ 
27r 


|l_ J^sin^ (^fcnr)(l-e-2«^'=)| . 
L 1 + Ofc J 


1 + ft 


( 12 ) 


This object (see also Eq. (HH) is the so-called cumulant generating function and 
coincides [U with the expression of Smacchia et al. [23]. We observe also that, in 
the limit s —>■ oo, we recover the expression for g^, the logarithm of the dynamical 
fidelity F{nT) = |('I'o|^'(?T-'r))|^ = discussed in Ref. |3D] (see Eq. (8) there). 


4. Steady state of the work probability distribution P^riW) for n ^ oo 


The question we are now going to address is if the probability distribution of the 
work Pt^fW) tends to “synchronize” with the periodic driving in the asymptotic limit 
tf —>■ oo. When viewing at the system stroboscopically at the times = nr (integer 
multiples of the period r), what we want to understand is if Pnr{W) converges towards 
a well defined asymptotic work distribution for n —> cx). The answer is positive, as 
we are now going to show. We know from previous work |28j that the quantum 


X The relevant quantity in Ref. m translates as follows in terms of our quantities; 


lft(nr)p 


2qk sin^ iukriT) 

[1 -I- Qk cos {2nknT)] 


We also mention that the analysis of Ref. |24| identifies the large-s limit of InG^^. ^{is)/L with a 
“surface” free-energy contribution 


C ^ dfc r "1 

-2/aurf = — In -t |j/fe(nT)PJ 


which coincides with our Qn- 









average of the work (i.e., the first moment of PnT{W)) indeed reaches a “steady 
state” for n —)■ oo, in the thermodynamic limit. [§] Moreover, as we have recently 
shown in Ref. [3^, the 5{W) Dirac-delta part of the distribution Pao{W) (see, for 
instance, Eq. (22) below) alias the large-s limit = lims_>.oo In (*s)/E — which 


corresponds to the dynamical fidelity — also reaches a well defined “steady state” for 
n —>■ oo: lim„_>oo Qn = 5oo- Our claim now is that all the cumulants of PnriW) reach 
such a “steady state”, and therefore so does the whole probability distribution. To 
see this, it is enough to show that the whole cumulant generating function, Eq. ( |l^ , 
reaches a steady state. With an argument which generalizes that of Ref. [3D], whose 
details are given in [Appendix A one can show that when n —>■ oo this quantity tends 
to the stationary value 


T=0 


InG^ °(*s) _ ,. 'i^G'nr 

- T - = - f 

h n^OO L 


{is) 


= 2 


In 

1 + \/l ~ f.k{s) 

Jo 27r 

1 

to 


(13) 


where 


a(s) = 4|r+| \r^\ . (14) 

The numerical results shown in Eig. [^perfectly confirm this analytic prediction: there 
is a transient — whose details depend on the parameters, for instance on the average 
field h(j — which then leads to an asymptotic result for n —>■ oo given by the simple 
analytic expression in Eq. 13 Fig.j^also illustrates (bottom panel) a case in which the 
frequency wq is such that Jo{2A/ujo) ~ 0, where Jo{x) is the Oth-order Bessel function, 
and there is coherent destruction of tunnelling [311133] : observe that has a 

very small magnitude for this value of ojq. 

The cumulants of the asymptotic work distribution are obtained as: 


Km = (-1 }"'-—lnG^='^(is) , 


(15) 


The first cumulant (which is the quantum average of the work performed) is given by: 


K^ = (W}^=4L - |r 


Pi -P f 

Pfc Ek ■ 


(16) 


This coincides with the result we would obtain by evaluating 

dfc 


lim / — 

n—>oo /„ 27r - 


{tpkiriT)] Hk{0) lifkinr)) - Ef (0) 


which is quite easy to calculate directly. The second cumulant is the variance of the 
work distribution and is given by: 


2 rdk r 


1 + 3 r 


Et 


(17) 


We notice that the P{W) tends to become narrower and narrower in the 
thermodynamic limit, as expected, because CToo/ (1+)oo 1/a/+- 

§ Strictly speaking, when L ^ oo the P{W) becomes narrower and narrower, on the scale of the 
average work {W) ~ L, with fluctuations which scale as y/L. Nevertheless, when L is large but not 
oo and the approximation of having set L ^ oo in transforming the sum into an integral holds only 
until a certain finite time t* ~ L, the question we are asking is meaningful, provided the “steady 
state” is effectively reached before t*. 




















Figure 1. Plot of ygj.gjjg g fQj. different n, compared to ^ for a 

Ising model with periodic driving h(t) = ho + ^cos (uiqI + <j>o]\ here hg = he = 1, 
A = 1, and two different values of the frequency are shown, too = 2 (top) and 
too = 0.3623 (bottom), corresponding to a situation in which Jo(2A/iuo) ~ 0 and 
there is coherent destruction of tunnelling (see discussion in the text). Numerical 
data are for L = 1000. 


5. Universal edge singularity at small W in Poo{W) 


Inspired by the results of Ref. [2l], we now discuss the behaviour of the asymptotic 
work probability distribution at small values of W, especially in connection with 
aspects which are independent of the details of the specific driving protocol. From a 
technical point of view, the small-IF behaviour of Poo (IF) is encoded in the large-s 
behaviour of Gao(is), which we can evaluate by means of Eq. (13|. We will show that, 


indeed, the important ingredients are: i) the value hi of the initial transverse field 
h{t = 0), and ii) the value of the average field hg (and the frequency Wq), determining 
if the Floquet spectrum shows a resonance, and is gapless, at fc = 0, or not. Indeed, the 
value of hi determines the position Wth of the singularity in P(IF) which we observe. 



















while the form of this singularity is determined by the possibility that the Floquet 
spectrum has a resonance at fc = 0, which is entirely determined by the time-averaged 
field ho (see Eq. (§) hitting a non-equilibrium quantum critical point h^i = hc-tlujo/^, 
with I integer, i.e., 


2\ho — hc\ = lojQ for some positive integer I 


(18) 


is obeyed. This conclusion holds for a generic periodic driving h{t): details can be 
found in Appendix C The small-VF universal behaviour of Pao{W) we describe below 
relies, in the end, only on the properties of the small-A: modes, in particular on the 
small-A: behaviour of \rt\'^ in Eq. (14): we find that if the resonance condition in 


Eq. (18) is fulfilled [Ji] then 


otherwise 


\rt? = \-\pk + 0{e)- 




(19) 


( 20 ) 


The precise values of a and jS depend on the specific form of driving, but the functional 
form of and the Floquet spectrum being gapless at fc = 0 or not depend only 


on the fulfillment of Eq. (18). Depending on fc; and ho, we can distinguish essentially 
three different behaviours of F’oo(bF) in the small-IT limit: 

Case a) hi he and non-resonant Floquet spectrum. In this case we obtain (see 
Appendix B for the derivation) an approximate analytical formula for Gao(is), valid 
when s S> l/|fci — hd- 


Goo(*s) ci 



-2s|/li-/le| 

g3/2 


( 21 ) 


The inverse Laplace transform predicts that the small-IF behaviour of Poo (IF) is given 

byH 


Poo(I^)-e^®“ 


^(IF) 


+ ^ ^/W-2\hi-he\9{W-2\hi - fcol)) , 


( 22 ) 


which applies whenever IF < (2 + ln2)|fci 
given by: 


a = 



— hc\. In both expressions the constant a is 


\hi-hdY^^ 

hi ) ’ 


(23) 


where is such that Ir^Tf — j't for small fc (see Eq. (20)). Eq. (22) predicts 
an edge singularity in the asymptotic work distribution function at a precise value of 
IF which is totally independent of the details of the periodic protocol (and even of 
the frequency) but depends only on the initial value hi of the field. The details of 


there are situations in which Eq. | |20| l is valid for and not for 

|r^p. Nevertheless, all the results are unchanged: since -f — t: exchanging the roles of 

|r^| and |r^| has no effect, as Eq. is symmetric under such an exchange. 

^ Remarkably, the form of the singularity in Eq. \22\ for the asymptotic work-distribution function is 
the same found in m in the case of a generic quench starting and ending into the same paramagnetic 
or ferromagnetic case. This is exactly what we are doing in this periodic driving protocol. 


As we show in 


Appendix C 













the protocol enter into the strength of the singularity (the coefficient a). We notice 
that the threshold 2\hi — is the energy which has to be provided to the system to 
generate an excitation in the k = 0 mode; moreover, also the form of the singularity 
is only determined, through the constant a, by the small-fc Floquet modes, as detailed 
in [Appendix C[ So, the behaviour at small W of the work distribution function is 
dominated by the modes of lowest energy: in retrospective, this is a very reasonable 
finding. We stress that the previous analytical expressions are approximations valid 
in a precise range of s or W. The only condition for the validity of these approximate 
formulas, as detailed in [Appendix B[ is that the driving field does not start from the 
critical point value, i.e., hi ^ he and there are no resonances at fc = 0 in the Floquet 

spectrum. We show some instances of the validity of Eq. (21) in the upper panel 

lnG„ 


of Fig. where we numerically evaluate for h{t) = ho + Acos{ujot + 4>o) 

/ versus s for different values of 


and plot R{s) = 
ho /ici hi ^ he 


In Goo (2s) 
L 


and wq. We see that, for large s, R{s) tends towards a constant. 


which Eq. (21) predicts to be a, Eq. (23), denoted by horizontal lines (the values of 


were obtained by a numerical fitting of for small fc, according to Eq. (20)). 


Overall, we see that there is a good agreement with the asymptotic value of R{s). 

Case b) hi ^ he and resonant Floquet spectrum (Eq. fulfilled). In this case 
the Gaussian analysis performed in [Appendix B[ fails. We find, nevertheless, an 
approximate form of Goo (is) for large s (s ^ l/\hi — he\) as 


Goo{is) ~ + Loc 


56 


— 2s\hi — hc\ 


■)' 


(24) 


as long as there is no coherent destruction of tunnelling, i.e., for a resonance of order 
I, we have that Ji(2A/uJo) ^ 0 {Ji being the Bessel function of the first kind of order 
I [H]). The resulting P{W), after inverse Laplace transform, has a very singular 
contribution 


Po 


,{W)c:ie^^°°(d{W)+aeLS'(W-2\hi-he\)0{W-2\hi-he\) + ---'^ . (25) 


The lower panel of Fig. shows plots of Re{s) = 


- goo] / [se-2|^.-ib] for 

different resonant cases mth I = 0 and ^ = 1. 

Case c) hi = he: When the initial field is critical, is gapless and the behaviour 
of Goo{is) is power-law rather than exponential. The upper panel of Fig. [^illustrates 
several cases with hi = he, for a driving of the form h(t) = ho + Acos(a;ot + 4>o) with 
ho such that the Floquet spectrum is resonant with I = 0, i.e., ho = hi = he (obtained 
for (fo = i7r/2), showing a clear power-law decay of the form 1/s: 


Goo (is) ^ (^l + 



(26) 


(Remarkably, this 1/s decay is valid also for uiq such that there is coherent 
destruction of tunnelling, for instance when wq = 0.3623..., where Jo{2A/u)o) = 0.) 
Correspondingly, we predict a step-singularity in Poo{W) 


Poo{W) ~ e^9“ (<5(W) + D e{W)) , 


(27) 


where 0 is the Heaviside function. We notice that in Ref. [21] the authors find a very 
similar formula for Pt{W) (when the time t is hnite) in an Ising chain undergoing a 




















Figure 2. Asymptotic cumulant generating function for a periodically driven 
uniform quantum Ising chain. In these figures we assume a driving h{t) = 
/i0 +Acos(a;o^ + 0o)? with A = 1 and (po = 0. (Top) Driving with ho ^ he = 1 and 
h-^ ^ he = 1, plot of R{s) = ! |■g-2|/li-l|syg3/2j por 

large s, we see a convergence towards a finite limit (a, as defined in Eq. ( |23| l and 
calculated by numerically fitting the coefficient oP" appearing in = cx^}p jA. 

for small k). (Bottom) Driv ing with hi ^ he = 1 but ho fulfilling the Floquet 
resonance condition Eq. ( |18[ ) (we take instances with resonances at / = 0 and 
I = 1). Plot of Rc{s) = versus s, showing the 

convergence towards a finite limit for large s, in agreement with Eq. | |24[ |. Notice 
that cases where there is coherent destruction of tunnelling fail to be described 
by such a formula. 





















Figure 3. Asymptotic cumulant generating function for a periodically driven 
uniform quantum Ising chain. In these figures we assume a driving h{t) = 
ho + Acos(cjo^ + 00)5 with A = 1. (Top) Instances of the critical case with 
h[ = ho = he (obtained taking 0o = and different frequencies: a power-law 

decay like 1/s is confirmed in agreement with Eq. ( |26| l. (Bottom) The critical 
case hi = he = 1 for ho = 0 and (po = 0. We see a 1/s^ decay with 6 = 3 or 1, 
depending on the Floquet spectrum being resonant (2 = lojo) or riot in k = 0. 


generic (non-necessarily periodic) driving which ends at the critical point. In our case, 
we consider the asymptotic behaviour and, thanks to periodicity, whenever hi = he 
the system not only ends but also starts at the critical point. It turns out that 
the case hi = he is quite rich, and other power-laws are possible if the Floquet 
spectrum is resonant at A: = 0 with non-zero values of the integer I in Eq. |T8| ). 
For instance, the bottom panel of Fig. shows cases with hi = he = A = 1 
and ho = 0: if the Floquet spectrum in A: = 0 is resonant (2 = lojo^ Eq. @) 
then decays like 1/s^, otherwise it decays like 1/s. A 1/s^ behaviour, 

Goo (is) ~ (1+ ^ j results into a mild quadratic increase at small W for the work 











distribution: P^{W) ~ (siW) + e{W)DW^l2j . 

6. Finite temperature results: irreversible entropy generation 

One can also calculate the average work performed in the finite temperature case, 
where the initial state is an equilibrium Gibbs state at a finite temperature ksT = j3~^ . 
Using Eq. for InGm-(M), and the fact that the average work performed in a time 
tf = nr is given by —idlnGuTiu)/du\u=o, one can readily arrive at the expression: 

ciA^ 12 I 12 

{W)nT = 'iLj — (1 - cos(2/ifcnr)) Irj!"! |r^| tanh(/3Ufc) , (28) 

which generalizes Eq. ( [I^ to hnite T and finite n. When n —)■ oo, the rapidly 
oscillating term cos (2^fcnT) gives a contribution that averages to zero, and we get: 

dk 2 2 

{W)oo=dLj — |r^| |r^| Ek ta.nh {13 E k) ■ (29) 

The statistical nature of the work for a hnite system — evidenced by Eq. 0 — calls 
for a second law of thermodynamics written as {W) > AF, which relates the average 
work performed to the difference in equilibrium free energy AF corresponding to 
the initial and hnal parameter values: here the equality holds only for a quasi-static 
process. Calling (W''''') = {W) — AF the difference between {W) and AF — and 
viewing it as the average irreversible work done —, one can recast the second law in 
the form (VU”'') > 0. The heat transfer between the closed system and the bath being 
zero, the entire contribution to the entropy generation is in fact due to (VE’'’’'), and 
one can dehne the irreversible entropy generated as A5'“''' = /3(IU)'“''' [THl[TH] . 

Since at stroboscopic times = nr the Hamiltonian returns to the original 
value H{0), so that AF = 0, one can dehne an irreversible entropy increase, in the 
limit n —)■ 00 , as AS'*’'’' = {3{W)ao- In Fig. 0 we show that for a driving protocol 
h{t) = 1-1- cos(wot), AS’"' indeed saturates to a steady state value, like the residual 
energy [35], displaying a sequence of well dehned dips and peaks: In the small loq 
regime, there are dips at certain frequencies for which Jq{2A/ujq) = 0, a consequence of 
the coherent destruction of tunnelling |41j . In the intermediate range, on the contrary, 
one hnds peaks at ujq = 4/p, with p integer, due to quasi-degeneracies in the Eloquet 
spectrum [35] . 

7. Conclusion 

In conclusion, we have studied a periodically driven transverse-held Ising model and 
analyzed the behaviour of the stroboscopic characteristic function Gnriis), and hence 
the stroboscopic work distribution function F„,-(tU), after n complete periods of 
driving. Our study establishes that, in the thermodynamic limit, Pnr indeed converges 
for n —>■ oo towards a well dehned steady state value Poo{W) which reproduces the 
exact asymptotic value of the hrst cumulant {W)co (i-O., the asymptotic value of the 
average work performed on the system) derived earlier |28j . In the limit s —^ oo, on 
the other hand, G„t(*s) reduces to the stroboscopic dynamical hdelity [50] . 

Eor large s, we are able to provide asymptotic analytical expressions for Goo (is) 
and, by means of inverse Laplace transforms, we can derive corresponding expressions 
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Figure 4. Figure showing that for a periodic driving protocol h{t) = 1 + cosojot, 
saturates to a steady state value in the asymptotic limit (n oo) like the 
residual energy with similar dip and peak structure in the low and intermediate 
a;o, respectively |28|. In the main figure peaks occur at loq = A/p while in the 
inset we show that dips in are located at Jo(2/a;o) = 0 for small loq. Here, 

^ = 10 and L = 1000. 

describing the small-TV behaviour of PaoiW). The small-W properties of Poo{W) 
depend strongly on the fact that there is a static critical point he, and any periodic 
driving induces further non-equilibrium critical points where the gap in the Floquet 
spectrum closes up. This finding is in line with the study reported in Ref. [33], where, 
however, the relevant regime was one of large-amplitude driving at large frequencies, 
and a rotating wave approximation was appropriate. Here, on the contrary, the 
exact resonances we find reside at low frequencies: for a fixed average held hg, at 
Wo = 2|ho — hc\/l. In any case, the form of the singularity in the work distribution 
turns out to be a useful detector of such non-equilibrium phase transitions. According 
to the way the external periodic driving held relates to these critical points we can 
observe different phenomena. The time-averaged value of the held and its initial 
value hi happen to be crucial. Whenever hi is different from the static critical point 
he and hp differs from any non-equilibrium critical point (i.e., 2|/io — hc\ ^ ^wq), the 
asymptotic F’oo(kF) is characterized by a universal edge singularity whose position 
Wth = 21 hi — he I depends only on the value of hi. The information about the specihe 
protocol appears only in the strength of the singularity. The low-energy behaviour of 
Poo{W) is essentially determined by the lowest excitation modes and the corresponding 
energy gap. If hi ^ he but hg is dynamically critical (i.e., 21hg — he| = Iloq for some 
integer 1), the form of this singularity changes, and becomes a Dirac delta derivative. 














A completely different behaviour emerges if the initial field is critical, h{ = he- in 
this case, the cumulant generating function (more precisely, hag a 

power-law decay, as 1/s, resulting in a step-function contribution Poo{W), but other 
decays, as 1/s^, can also be seen if the Floquet spectrum is resonant with I ^ 0. 
Generalizing our investigation to the finite temperature case, we have shown that the 
irreversible entropy AS'*’’*', obtained using the first cumulant of the finite temperature 
characteristic function, also synchronizes with the periodic driving for n —> oo and 
converges to a steady state value for large ojq- 

Summarizing, we see a strong relationship between the features of Poo(hF), the 
existence of a critical point and the way the driving field relates to it. This work is a 
first step towards the application of time-periodic probes to understand the existence 
of a quantum phase transition by looking at the work distribution function. In this 
sense, we are generalizing the very interesting works in Refs. [2Ql EH EH |43], which 
refer to the case of a sudden quench. In this perspective, it is also interesting to see 
if it is possible to induce non-equilibrium phase transitions in systems without static 
transitions and how this influences the work statistics. Another possible direction 
is to see how the quantum driven system being regular or ergodic (m EH sa szi) 
influences the work statistics. A lot of work still remains to be done, starting from 
the consideration of other tractable cases, like the Dicke model [43] . 


Appendix A. 

In this appendix we show how the stroboscopic cumulant generating function tends 
towards the stationary value given by Eq. (13). Our first step is to expand the 

e fo] 


logarithm in Eq. (12), which we report here for the reader’s convenience 
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2tt 
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1 + <7fc 


sin^ (^fcnT)(l - e . (A.l) 


To that purpose, we first show that the second term inside the logarithm is < 1. We 
know that sin^(n/ifcT) < 1 and |1 — < 1 whenever s > 0 and s ^ oo. As for 

the overall prefactor, we notice that 


Cfc = 


1 + 


= 4 \r1 




(A.2) 


where the value = 1 is obtained for = 1/2. Hence, the expansion of the 

logarithm is certainly possible for all s < oo. Defining 


ik{s) = 


‘^Qk 


1 + 


(1 - (1 - 


we have: 
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sin (^fcnr) , 


(A.3) 


(A.4) 


where we have exchanged the integral and the sum over m, due to the dominated 
convergence theorem. We then write a binomial expansion of the sine term in terms 
of exponentials: 


sin^’” inkUT) = 


i-iy 


■E 

3=0 


2m 
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(-l)^e 




(A.5) 












The sum over j, for each m, has a finite number of terms: there is no problem in 
exchanging the integral over k with this sum. We now observe that the j ^ m terms 
contain rapidly oscillating factors and vanish in the limit n —> oo, thanks to the 
Riemann-Lebesgue lemma and the smoothness of the factors Hence, in the limit 
n —>■ 00 , we retain only the j = m terms and write 


lim 

n—>-oo 




T" V-Lp-IST). (A. 6 ) 
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We can write this expression in a closed form, by defining 

4 m 'yrn / m 
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and noticing that: 
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which can be integrated to give 
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(A.7) 


(A.8) 


(In the last steps we have substituted ^ = sin^(? 7 ) and rj' = rj/2.) In conclusion, we 
arrive at the desired result: 


lim 

n—>oo 


lnGj=0(7s) 


In 

1 + \/l ~ C/c(s) 

/o Stt 
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csi 

_1 


(A.9) 


where ^fc(s) is given in Eq. (A.3|. 


Appendix B. 


We analyze here the asymptotic large-s behaviour of We start by assuming 

that the initial transverse field is not critical, /li 7 ^ /ic = 1 , so that there is a gap in the 
spectrum of the initial Hamiltonian, Ek > |hi — I| > 0. We need this condition 
because we would like to expand Eq. (A.9) to lowest order in e“®^'' and this is 


possible, provided Ek > 0 and s S> I/|/ii — 1|. With the previously defined shorthand 
^k = ^k{s —>■ 00 ) = 2qk/{l + qk) (which is a positive quantity in [0,1]) we can rewrite 
Eq. (A.9) as 


lnG^=0(*s) ^ 

In 

1 + Vl - 6 + 

L “j 

0 2. 

2 


(B.l) 


Expanding this to first order in e (with the assumption Cfc 1, that is 

P 7 ^ 1/2, see comments below) we find: 


In G^-°(js) ^ logGoo(is-)> 00 ) 
L ~ L 
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^k 




2 ^ (1 + 
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The expansion is questionable whenever there are fc-points such that = 1. We 
will see below that this is indeed the case at fc = 0 whenever the field oscillates 


around a non-equilibrium critical value given by Eq. (18) and the Floquet spectrum 
is resonant in fc = 0. Even restricting ourselves to non-resonant cases, we would 
possibly find fc-points where ^ko = but this time with ko > 0. This would still 
seem to be an issue at first glance: indeed, near these points we would expand 
quadratically, = 1 — A^(fc — fco)^ + 0{{k — ko)"^), and we would therefore have, 
in the integrand of Eq. (B.2), some logarithmic singularities originated by terms of 


the form 1/ {k — ko). In the neighborhood of these points where ^ko = 1? however, a 
different expansion is more appropriate. Considering one such ko where ^ko = !> we 
can write the contribution to lnG'^'^{is)/L from the neighborhood of such a point, 
see Eq. (B.l), as 


'^o+^dk , 

— In 
27r 


'ko-e 


l-he 


-sEk 


2 In 2 


^_^Q-2sEkg f ^^-2svkgk' 

./_« 27r 


where we have approximated ^k — 1, expanded the logarithm and introduced the group 
velocity Vk^ = dEk/dk\f,^. This is indeed a convergent expression. More importantly, 
when ko ^ 0 the singularities needing such a special treatment are in a region where 
the integrand is exponentially smaller, due to the prefactor then the main 

contribution coming from the fc = 0 region, which we are now going to analyze. 

To proceed with the case fci ^ fc-c = 1, we expand Ek up to second order in fc 


Ek — \hi — 1\ + 


hi 

2\hi-l\ 


0{k‘^) 


(B.3) 


Because we are assuming |/ii — 1| > 0 and s S> l/|fci — 1|, only the fc such that 
fc ^ \/\hi — l|/fci contribute significantly to the integral in Eq. (B.2|: we would 
be willing, therefore, to approximate the integral with a Gaussian one, extending the 
upper integration limit to oo and expanding the factor multiplying to the lowest 

order in fc. But here the fc = 0 point plays a tricky role. For a generic periodic driving 
where h{t) oscillates around the average value ho one quickly realizes, by focusing on 


the evolution operator Ufe(r, 0) for modes with a small fc (as detailed in Appendix C) 
that there are two cases: 


i) The Floquet spectrum is resonant in fc = 0 whenever 

2\ho — hc\ = liOo for some positive integer I 


and then fj,^ (x ±k and 


.+ |2 


= l-^4+oiP) 


(B.4) 


(B.5) 


hence = 1 ^^od the expansion in Eq. (B.2) is inappropriate. As detailed in 
[Appendix C for a sinusoidal driving of the form h{t) = ho + Acos{ujot + (l)o) 
there are special situations where there is coherent destruction of tunnelling 
(CDT) imilT]: if Ji{2A/uJo) = 0, for the value of I that realizes the resonance, 
then a gap is opened in the Floquet spectrum and = 0; posing no problem 


with Eq. (B.2). 


















ii) If, on the contrary, the Floquet spectrum is not resonant in fc = 0, then 


\rt? = 


a 


o{k^). 


(B.6) 


hence = Oj again posing no problem with Eq. (B.2). 


So, restricting our consideration to non-critical initial helds, /i, ^ he, and 
non-resonant Floquet spectrum, and performing the appropriate Gaussian integral 


emerging from Eq. (B.2), we finally arrive at: 

T=0/-„\ „,2 


InG, 
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where we have defined 


16 y^ V h[S 

while a is the 
. Since s is large, we can equivalently recast this 

(B.8) 

3/2 


where we assumed s ;i> l/|hi — 1|, while a is the constant appearing in the quadratic 
expansion of \r^\^, see Eq. 
equation as 

aL 


1+ -h-' 


a 


a = 




16 v^ \ K J ■ 

Performing the inverse Laplace transform, it is not difficult to show that the P(W) 
associated to Eq. (B.81 is: 

Poo(VF) + ‘^^w-2\hi-l\ 0{W-2\hi-l\) + ---^ . (B.IO) 

So, our theory predicts a square-root edge singularity in the asymptotic work 
distribution at a precise value Wth = 2 |/ii — 1| of W, i.e., simply the initial gap of 
the system, a value which is totally independent of the details of the periodic protocol 
(and even of the frequency), which only enter into the prefactor a. 

Finally, let us briefly consider the remaining cases where the previous Gaussian 
analysis fails. There are, essentially, two cases left: 

(i) Case ^ he and resonant Floquet spectrum. The starting point of this case is 
Eq. (B.l), with given by Eq. (B.3) but with the complication that ^k-^o = 1 
in the vicinity of fc = 0. An analysis of the singularity emerging shows that the 
leading term for large s is now given by: 


lnG^=0(*s) 


- 




Exponentiating, we can equivalently rewrite: 

G^=°(is) ~ e^3“ (l -k ae Ls -p • • •) . 


(B.ll) 


(B.12) 


The resulting P(W), after inverse Laplace transform, has a very singular 
contribution proportional to a Dirac’s delta derivative: 

Poo(VF) ~e^3”y(IF) -f UeL S'(W - 2\hi - 1\) e (W - 2\hi - 1\) -k---) 

(B.13) 

















(ii) Cases with hi = he- Whenever the initial field hi coincides with the critical field 
he the initial spectrum is gapless, and this changes completely the large-s 
behaviour of G'^^{is), from exponential to power-law. The scenario is quite 
rich. For instance, when ho = he the behaviour is of the form 


leading to a small-VF probability distribution 

Poo(VF) (^(VF) + D 0(W)) , 


(B.14) 


(B.15) 


where 0 is the Heaviside function. But this does not exhaust all the possibilities: 
when ho = 0 we find that the leading asymptotics of Geo (is) is 1/s^ rather than 
1/s whenever the Floquet-resonant condition Eq. (18) is fulfilled. A thorough 
study of the gapless scenario is left to future studies. 


Appendix C. 


In this appendix we prove the statements leading to the resonance condition Eq. (18) 


and the related small-fc expansions Eqs. (B.5) and (B.6). Let us start with the 
resonance condition. The argument is very similar to the one reported in the 
Supplemental Material of Ref. |25j. Let us consider Eq. ([^ with a generic time- 
periodic h{t) with period r = 2'k/ujo and apply to it the time-dependent rotation 
Vfc(t) = exp {—if{t)az), where f{t) = dt' {h{t') — ho). The matrix Mk(t) in the new 
representation has the form 


hit) = Vl(t)Hfc(t)Vfc(t) - iWlit)Ykit) = 


ho — cos k —i sin k 
ismke~^^t'{t) _/jp_|_eosfc 


(C.l) 


The unitary matrix Vfe(t) has the remarkable property that Ykir) = Vfe(O) = 1. 
Because of that, we can immediately write the time-evolution operator over one period 
T in the original frame as 


Ufc(T,0) = 


(C.2) 


In essence, the Hamiltonian in the rotated frame directly determines the form of 
the Floquet modes and quasi-energies, obtained by diagonalizing l[Jfc(T, 0). With 
a vanishing driving, the Floquet spectrum is given by the eigenenergies of the 
Hamiltonian folded in the first Brillouin zone [55], therefore we have a resonance 
whenever 2Ek = lojo for some integer 1. For a finite external field hit), most of the 
resonances turn into avoided crossings, but for the modes with fc = 0 and k = tt. 


Indeed, we can see from Eq. (C.l) that the field couples to sinfc, which vanishes at 


fc = 0,7r. If we diagonalize Eq. (C.2) at those values of k, we easily find the single¬ 
particle Floquet quasi-energies as = ±|/io — 1| and /rj = ±|/io -I-1|. The resonances 
at fc = 0 are particularly important since, for s —> oo only the smallest values of Ej. 


matter (see Eqs. (13) and (HI)), and these occur near fc = 0. Folding the quasi-energy 
into the first Brillouin zone, we can see that the Floquet spectrum is resonant at fc = 0 
when: 


2|fco ~ 1| = 1^0 for some positive integer I . 


(C.3) 
















We now move to establishing Eqs. (B.5) and 


to consider that the overlap factors |r? 


For this purpose we have 
= K'/'fc (0)|V^f®)P are obtained from the 


Floquet modes |<^*(0)), which are the eigenstates of an Hermitian operator, the 
Floquet Hamiltonian dehned as 


g-irHf = dtHfe(t) ^ 


(C.4) 


(The Floquet quasi-energies are the eigenvalues of H^.) From the above discussion, 
see Fqs. (C.l) and (C.2), we immediately see that 


Co 


ho-I 0 

0 1-ho 


Because the Floquet quasi-energies are dehned up to translations of an integer number 
of Wo, this Hamiltonian is equivalent to 


Co 


h 

0 —h 


(C.5) 


where h = (/iq — 1) —Zwo/2 is /ig —1 folded in the hrst Brihouin zone [— wo/2, wo/2]. We 


see that h — 0 whenever there is a resonance, i.e., Eq. is valid. Eq. (C.4) shows 


that IHI^C must have the same symmetry properties of IHIfc(t). We see, from Eq. (7), 
that the symmetry relation 

= (T^Wk(t)a^ 

is valid at ah times. Because az is a time-independent unitary transformation. 


Eq. (C.4| implies that = az^^ Uz- Because of Eq. (C.5), and the relations 


Tr[IHIfe(t)] = 0, = az, azaxOz = —ax, azayaz = —o'y, we can write the following 


second-order-in-fc expansion □ of 


■‘'■small k 
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In general the coefficients ax, ay and Uz are non-vanishing; they can vanish in some 
cases, giving rise to interesting phenomena which we will discuss later. Whenever 


the resonance condition Eq. (C.3) is not fulHlled (hence h ^ 0), the second-order-in-fc 
expansion of Floquet modes (expressed in the basis H = {|0), cj'c^^ |0)}) is 
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which applies to the case h > 0; these two states should be exchanged if h < 0. 
Moving now to the initial Hamiltonian ground state It/'®’*), the diagonalization of 


■*" The vanishing of the trace of ]H[(( at any k comes from the vanishing of the trace of Hj. (t) and the 
formula Tr[El((] = F JJ" Tr[Hj, (t)]dT, which is a corollary of the Liouville’s theorem m- 



















Eq. 0 immediately gives (for hi> \) 


[IV'fmall k)\B ~ 


Hence, for the overlap we find 
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which is indeed Eq. (B.6); this formula is valid also in the case h < 0 and hi < 1. 
If h < 0 and hj > 1, or h > 0 and hi < 1, we find + 0{k^), but the 


(see Eq. (A.2)). 


crucial ingredient determining Eq. (A.9) is identical, since — a^k"^ in both cases 
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For the resonant case h = 0 we find 
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which (if hi > 1) gives rise to 
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We notice that this formula is 


V 


This is indeed Eq. (B.5) with /3 = a^/ + 

valid if Ux — iay ^ 0. 
coherent destruction of tunnelling (CDT) [55]. where = a-y = 0, and we fall 
back to Eq. (B.6). In the Supplemental material of Ref. [55] we show that, in the 


y -f- Kj. This is generally true, up to special cases where there is 


case of a sinusoidal driving h{t) = hg + Acos{ujot + 4>o), CDT occurs if hg = 1 and 
Jo(2A/a;g) = 0. More in general, if the resonance condition Eq . ([C.3[ ) is valid for I ^ 0, 
we can show, exactly with the same arguments used in Ref. |28j. that there is CDT 
whenever Ji{2A/ujq) = 0. 
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